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Abstract 


A  new  method  for  characterization  of  unknown  targets  using  passive  multispectral 
polarimetric  imagery  is  presented.  Previous  work  makes  use  of  a  pBRDF  derived  equation 
for  the  degree  of  linear  polarization  and  with  the  aid  of  measurements  at  multiple  incident 
angles  estimates  refractive  index  and  reflection  angle.  This  work  uses  known  incident 
and  reflection  angles  along  with  dispersion  equations  and  polarimetric  data  at  multiple 
wavelengths  to  recover  the  index  of  refraction.  Although  imagery  is  collected  with  a 
division  of  time  polarimeter  and  a  spectral  filter  wheel  in  iterative,  manual  steps,  the 
new  algorithm  could  be  applied  to  any  set  of  registered  multispectral  polarimetric  images 
most  notably  those  produced  by  a  recently  introduced  division  of  focal  plane  multispectral 
polarimetric  sensor.  Experimental  results  are  presented  showing  the  novel  algorithm’s 
ability  to  classify  and  characterize  a  range  of  materials. 
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MATERIAL  CHARACTERIZATION  USING  PASSIVE  MULTISPECTRAL 


POLARIMETRIC  IMAGERY 

I.  Introduction 

Optics,  the  study  of  light  and  its  interaction  with  matter,  has  a  rich  history  dating 
back  to  the  first  use  of  mirrors  many  centuries  before  the  common  era  [12] 
and  obligates  the  mention  of  great  scientists  and  mathematicians  such  as  Snell,  Fermat, 
Huygens,  Newton,  Faraday,  Maxwell,  Planck,  and  Einstein  [16].  The  theories,  experiments, 
and  equations  generated  by  these  great  minds  (and  many  others)  allow  our  modern 
understanding  of  light’s  duality;  its  existence  as  quantized  particles  called  photons  and  as 
waves.  Of  particular  interest  to  this  research  is  the  propagation  of  light  shown  by  Maxwell’s 
equations  to  take  the  form  of  a  Transverse  Electromagnetic  (TEM)  wave. 

Polarization  (discussed  further  and  graphically  illustrated  in  Chapter  II)  is  defined  as 
the  orientation  of  the  electric  field  vector  of  a  TEM  wave  with  respect  to  its  propagation 
axis.  The  polarization  state  of  light,  although  imperceivable  to  the  human  eye,  can  be 
mathematically  derived  and  used  to  determine  valuable  information  about  the  object  from 
which  the  light  is  scattered.  The  goal  of  this  research  is  to  use  polarization  data  at  multiple 
wavelengths  in  a  novel  way  such  that  an  unknown  object’s  index  of  refraction  is  recovered. 

Over  the  past  three  decades  efforts  have  increased  to  “fuse”  data  from  sensors  of 
differing  modalities  (Multispectral  (MS),  thermal,  etc.)  to  more  accurately  identify  targets 
of  interest  [50].  Advances  in  nanofabrication  now  allow  construction  of  a  unified  MS 
polarimetric  sensor  on  one  substrate  [19],  eliminating  the  need  to  fuse  data  from  separate 
MS  and  Polarimetric  Imaging  (PI)  sensors.  The  algorithm  developed  in  this  work  could 
exploit  the  data  generated  by  this  new  sensor  to  more  accurately  detect  high  value  targets. 
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1.1  Problem  Statement  and  Hypothesis 

The  refractive  index  of  a  solid  material  is  a  measure  of  the  speed  at  which  light 
travels  through  that  material  relative  to  the  speed  of  light  in  a  vacuum  and  is  dependent 
upon  the  incident  light’s  wavelength.  Depending  on  the  Remote  Sensing  (RS)  technique 
employed,  recovering  an  object’s  index  of  refraction  may  require  the  determination  of  a 
large  set  of  unknown  and  often  coupled  parameters.  With  regards  to  polarimetric  RS  there 
exists  an  equation  (given  and  further  defined  in  Chapter  II)  that  expresses  the  Degree  of 
Polarization  (DOP)  of  light  scattered  from  an  object  as  a  function  of  the  object’s  refractive 
index  and  the  light’s  angles  of  incidence  and  reflection  [18,29,34,37,42],  Thilak  et  al.  used 
this  equation  to  estimate  an  object’s  complex  refractive  index  and  the  angle  of  reflection 
by  acquiring  multiple  images,  each  with  a  different  angle  of  incidence  [42],  This  research 
seeks  an  answer  to  the  following  problem  statement:  could  an  estimate  be  made  of  a  target’s 
index  of  refraction  with  a  single  MS  polarimetric  image? 

It  is  reasonable  to  assume  the  Sun’s  and  a  remote  sensing  platform’s  zenith  angle 
with  respect  to  a  target  (0,  and  6r  respectively  in  Fig.  1)  are  known  via  Global  Positioning 
System  (GPS)  and  inertial  measurement  unit  ephemeris  data  and  a  data  collection  time 
stamp.  Therefore,  a  novel  hypothesis  is  formed  that  with  known  0;  and  6r  angles  and 
DOP  measurements  collected  simultaneously  at  multiple  wavelengths,  the  target’s  index  of 
refraction  can  be  estimated  via  a  non-linear  least  squares  fit  using  two  dispersion  equations. 


Figure  1:  Viewing  Geometry  [18] 
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1.2  Scope  and  Assumptions 

In  order  to  remotely  recover  one  of  a  target’s  many  potential  physical  properties 
such  as  its  refractive  index,  a  method  must  be  found  to  isolate  that  parameter  in  a  set  of 
equations  such  that  solving  for  the  parameter  of  interest  is  feasible.  To  simplify  this  multi¬ 
dimensional  problem  assumptions  are  frequently  made  about  the  target’s  environment,  its 
other  physical  properties,  and  the  RS  viewing  geometry.  The  assumptions  necessary  to 
ensure  successful  completion  of  this  work,  thus  limiting  its  scope,  are  detailed  below.  Many 
of  these  assumptions  align  with  those  made  by  Hyde  [18]. 

The  first  two  assumptions  impose  limitations  on  the  target’s  environment.  The 
illumination  source  is  assumed  to  be  passive  and  unpolarized  (or  randomly  polarized)  and 
the  target  is  expected  to  be  stationary  with  respect  to  its  environment.  In  an  operational 
environment  the  first  assumption  simply  means  the  Sun  is  expected  to  be  the  light 
source.  This  condition  allows  for  a  common  reduction  in  scope  related  to  calculating  the 
polarimetric  state  of  the  light  scattered  off  the  target  as  described  in  sections  2.2.3  and  2.3.3. 
A  stationary  target  is  important  as  pixel  misregistration  can  cause  misclassification  and/or 
mischaracterization  of  the  target.  Two  components  of  this  work’s  experimental  instrument, 
a  manual  spectral  filter  wheel  and  a  focusing  lens  in  front  of  the  camera,  have  the  potential 
to  contribute  to  pixel  registration  issues.  Part  of  this  work’s  data  analysis  algorithm  includes 
averaging  experimental  values  on  a  pixel-by-pixel  basis  across  measurements  made  at 
different  wavelengths,  so  any  change  in  image  alignment  due  to  rotation  of  the  filter  wheel 
or  defocus  due  to  chromatic  aberration  apparent  on  the  camera’s  focal  plane  between 
measurements  could  cause  pixel  misclassification  and/or  mischaracterization.  The  method 
used  to  mitigate  this  issue  is  described  in  sections  3.1.3  and  4.2. 

The  third,  fourth,  and  fifth  assumptions  limit  the  physical  properties  of  the  targets 
of  interest  (dielectrics  and  metals).  All  targets  are  assumed  to  be  non-magnetic,  of  room 
temperature  (300K),  and  isotropic.  The  refractive  index  of  a  material  r\  is  defined  by  its 
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intrinsic  properties:  permeability,  conductance,  and  permittivity  [18,28].  Any  changes  to 
these  properties  must  be  accounted  for  in  assumptions  regarding  their  applicability  to  this 
work  or  by  the  functional  form  used  for  refractive  index.  The  materials  examined  in  this 
work  do  not  have  a  magnetic  response  at  optical  frequencies,  therefore  their  permeability  is 
assumed  to  not  be  of  concern.  Changing  conductance  does  not  apply  to  dielectrics  since  by 
definition  dielectrics  are  non-conductive.  The  conductance  of  metals  changes  as  a  function 
of  temperature  which  is  accounted  for  in  the  damping  term  of  the  Drude  model-derived 
functional  form  of  the  refractive  index  (see  Eq.  (3))  [16].  A  material’s  permittivity  (its 
susceptibility  to  dipole  moment  polarization  upon  interaction  with  an  external  electric  field) 
relative  to  that  of  free  space  is  its  relative  permittivity.  Given  the  assumptions  already 
made  here,  the  index  of  refraction  is  defined  as  simply  the  square  root  of  a  medium’s 
relative  permittivity.  The  results  of  this  work  display  changes  in  relative  permittivity  over 
wavelength  and  an  assumption  regarding  temperature  is  made  to  account  for  changes  in 
permittivity  due  to  temperature.  Most  available  refractive  index  values  are  experimentally 
found  at  room  temperature  [28]  and  the  targets  of  interest  exclude  those  of  high  temperature 
such  as  engine  exhaust,  so  the  assumed  temperature  of  300K  is  reasonable  and  allows  for 
comparison  of  this  work’s  results  to  those  available  in  the  literature.  In  general,  permittivity 
changes  over  the  surface  of  a  material,  but  even  at  visible  wavelengths  the  contributions 
of  a  large  number  of  atoms  are  averaged  during  imaging  and  represent  the  behavior  of 
an  isotropic  dielectric  [16].  Assuming  isotropy  allows  for  the  reduction  of  the  generic 
dielectric  tensor  to  a  scalar,  greatly  simplifying  the  mathematical  framework  necessary  for 
finding  the  index  of  refraction. 

The  final  two  assumptions  limit  the  viewing  geometry  and  imaging  instrument 
conditions.  All  imaging  is  assumed  to  occur  in  the  specular  plane  without  any 
pixel  saturation.  The  specular  plane  is  defined  as  that  containing  the  target’s  surface 
normal  vector  and  the  incident  and  reflected  light’s  propagation  direction  vectors.  The 


4 


experimental  instrument  used  in  this  work  is  currently  set  up  for  data  collection  in 
only  the  specular  plane.  Extension  to  out-of-plane  measurements  is  left  to  future  work. 
Pixel  saturation  causes  inaccuracy  of  the  calculated  polarimetric  data  and  leads  to  pixel 
misclassification  and/or  mischaracterization,  so  care  was  taken  to  ensure  no  pixel  saturation 
during  the  acquisition  of  irradiance  data.  Implementation  of  this  algorithm  should  be 
limited  to  systems  that  either  have  a  variable  source  intensity  or  a  camera  with  a  high 
enough  dynamic  range  to  ensure  no  or  very  little  pixel  saturation. 

1.3  Thesis  Overview 

An  introduction  to  the  problem  addressed  by  this  research,  the  author’s  proposed 
answer  to  that  problem,  and  assumptions  necessary  to  simplify  the  problem  were  presented 
in  Chapter  I.  Chapter  II  steps  through  the  academic  and  mathematical  background 
necessary  to  understanding  the  author’s  approach  to  proving  the  given  hypothesis.  Chapter 
III  shows  how,  given  the  background  developed  in  Chapter  II,  the  author  arrived  at  the 
stated  hypothesis  and  how  that  hypothesis  was  tested  via  experimentation.  Experimental 
results  and  analyses  thereof  are  given  in  Chapter  IV  while  Chapter  V  offers  concluding 
remarks  and  potential  future  extensions  to  this  work. 
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II.  Background  and  Related  Research 


An  attempt  is  made  in  this  chapter  to  ease  the  reader  into  the  presented  research  by 
starting  with  a  broad  concept  and  continually  narrowing  the  topic  focus  until  the  specific 
mathematical  framework  used  in  this  work  is  presented  at  the  end.  An  overview  is  given 
of  conventional  remote  sensing  approaches  and  the  case  is  made  for  the  specific  concept  of 
combining  multispectral  and  polarimetric  remote  sensing.  The  scope  is  then  focused  further 
by  providing  the  theory  underlying  each  term  of  this  work’s  title.  Finally,  the  chapter  delves 
into  the  mathematical  framework  developed  in  related  research  and  leveraged  in  this  work. 

2.1  Remote  Sensing 

The  term  remote  sensing  is  attributed  to  any  effort  to  obtain  information  about  an 
object  without  physically  coming  into  contact  with  said  object.  A  few  of  the  many 
possible  examples  of  remote  sensing  include  gathering  information  about  the  cosmos  with 
a  telescope,  mapping  Earth’s  weather  patterns  using  a  satellite,  or  measuring  the  speed 
of  a  moving  vehicle  with  a  radar  gun.  Each  of  these  examples  involves  the  detection 
and  measurement  of  Electromagnetic  (EM)  radiation  of  different  wavelengths  reflected 
or  emitted  from  the  object(s)  of  interest.  The  specific  wavelength  boundaries  of  each 
named  region  of  the  EM  spectrum  vary  depending  on  the  source  consulted  with  one  such 
interpretation  offered  in  Fig.  2  [38], 

wavelength  ,  ,  ,  ,  1  n  ,  1  cm  1  mm  ,  1  mil  |  1,|i  ,  ,  1  nm  1 A  |  1  |>m 

A(m)  103  102  10  1  10'1  10-2  10-3  IO*4  10-5  10-6  10-7  IO*8  10'9  10*10  io-11  io-12 


Figure  2:  Electromagnetic  Spectrum  [38] 
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Remote  sensing  of  EM  energy  can  occur  in  a  lab  setting,  but  most  often  RS  refers 
to  gathering  information  about  a  distant  object  from  Earth’s  surface  or  information  about 
Earth  from  an  aerial  or  space  asset.  In  order  to  detect  radiation  of  a  particular  wavelength, 
a  system’s  sensor  must  be  designed  such  that  it  is  sensitive  to  that  wavelength  and  can 
“see”  that  wavelength  through  Earth’s  atmosphere.  Energy  from  the  Sun  in  the  form  of 
mainly  Infrared  (IR)  and  visible  light  impinges  on  Earth’s  atmosphere  (Fig.  3a)  where  it 
is  absorbed  or  transmitted  according  to  Fig.  3b  [21].  The  wavelengths  at  which  energy 
is  blocked  correspond  to  specific  atmospheric  molecules  that  absorb  that  wavelength  of 
incoming  radiation.  Acquiring  RS  data  is  limited  to  the  non-blocked  spectral  regions,  the 
“atmospheric  windows.”  Figure  3c  shows  which  RS  systems  are  useful  in  each  atmospheric 
window  region.  The  next  section  gives  a  cursory  explanation  and  examples  of  but  a  few 
conventional  RS  systems  [21]. 


c: 

o 


— h  h-  Human  eye 

Photography  Thermal  scanners 
H - H  H - H 

Multispectral  scanners  Radar  and  passive  microwave 

H - H  H - H 

— 1 — i — 1  i  1 1 1 1 1 - 1 — '  '  i  i  m  1 1 - 1 — I  -  i  i  1 1 1 1 - — | - 1 - 1 - 1 - r 

0.3  /im  1 /i,m  10ju,m  100  /im  1mm  1m 
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(c) 

Figure  3:  Spectral  Characteristics  of  (a)  the  Sun  and  Earth,  (b)  Atmospheric  Transmittance, 
and  (c)  Common  Remote  Sensing  Systems  [21] 
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2.1.1  Conventional  Approaches. 

As  shown  by  Fig.  3c,  RS  systems  can  be  classified  by  the  wavelength  region  they 
operate  over,  but  of  greater  interest  to  this  work  is  the  classification  of  a  RS  system  as 
active  or  passive.  Active  systems  produce  a  signal  of  a  particular  EM  energy,  send  that 
signal  towards  an  object  or  scene  of  interest,  and  then  use  the  measurement  of  the  returned 
signal  to  make  observations  about  the  object.  Passive  systems  detect  naturally  occurring 
EM  energy,  such  as  that  from  the  Sun,  scattered  and/or  emitted  from  an  object. 

Examples  of  active  remote  sensing  systems  include  Radio  Detection  and  Ranging 
(RADAR),  Light  Detection  and  Ranging  (LiDAR),  and  Sound  Navigation  and  Ranging 
(SONAR).  RADAR  systems  transmit  radio  frequency  EM  waves,  often  times  with  a 
specifically  generated  polarization  state,  toward  a  region  of  interest  and  receive  and 
detect  these  EM  waves  when  reflected  from  objects  in  that  region.  Among  many  other 
applications,  RADAR  is  used  for  detection  and  tracking  of  targets  of  interest,  creating  two 
and  three  dimensional  maps,  collision  avoidance,  and  Earth  resources  monitoring  [31]. 
RADAR’s  advantages  include  all-weather  performance,  day  and  night  capability,  and  an 
ability  to  penetrate  a  modest  depth  into  the  surface  of  the  Earth.  The  biggest  disadvantage  of 
RADAR  is  a  system’s  resolution  (its  ability  to  discern  between  two  closely  spaced  objects) 
is  limited  by  its  antenna  size  and  wavelength.  A  technique  called  Synthetic  Aperture 
Radar  (SAR)  is  used  to  overcome  this  limitation,  but  for  the  sake  of  brevity  will  not  be 
covered  here.  An  example  of  a  RADAR  product  collected  over  two  radio  wavelengths  with 
multiple  polarization  states  displayed  in  false  colors  is  shown  in  Fig.  4a  [32], 

LiDAR  systems  emit  laser  light  in  the  near  IR  and  visible  portions  of  the  EM  spectrum 
and  measure  the  time  it  takes  for  the  sensor  to  detect  the  reflected  light.  The  laser  light  is 
emitted  from  the  system  in  rapid  succession  to  form  a  “point  cloud”  of  raw  LiDAR  data 
in  which  each  point  represents  a  three  dimensional  coordinate  with  respect  to  the  known 
sensor  location.  Return  intensity,  number  of  multiple  returns  (when  the  light  reflects  off 
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more  than  one  object  before  returning  to  the  sensor),  and  scan  angle  information  is  also 
encoded  as  part  of  each  data  point.  Some  applications  of  LiDAR  include  water  depth 
analysis,  terrain  mapping,  and  tree  canopy  profiling  [25].  One  distinct  advantage  of 
LiDAR  is  its  raw  data  is  georeferenced  from  inception,  simplifying  image  production  [21]. 
Calibrating  the  sensor’s  location  to  produce  the  georeferenced  raw  data,  however,  can  be 
a  tedious  process  since  a  surveyed  ground-based  GPS  location  near  the  survey  area  is 
necessary.  An  example  of  a  LiDAR  product  is  shown  in  the  upper  portion  of  Fig.  4b 
showing  the  “bare  ground”  elevation  of  Mount  Hood  in  Oregon  [27]. 


(a)  RADAR  -  Mount  Everest  [32]  (b)  LiDAR  -  Mount  Hood  [27]  (c)  SONAR  -  Lau  Basin  [23] 

Figure  4:  Active  Remote  Sensing  Product  Examples  and  Their  Optical  Counterparts 


SONAR  systems  rely  on  mechanical  acoustic  waves  for  navigation  and  object 
detection  purposes  rather  than  on  the  EM  waves  used  by  other  RS  systems.  Sound 
waves  travel  longitudinally  (versus  transversely)  in  most  mediums  via  compression  and 
rarefaction  of  the  medium  ahead  of  the  wave  in  the  direction  of  propagation.  The  sound 
wave  pulses  created  by  a  SONAR  system’s  transducer  reflect  off  objects  of  differing 
mediums  and  that  reflected  echo  is  received  by  the  system.  High  frequency  SONAR 
is  used  for  mapping  sea  floors,  medium  frequency  military  systems  detect  submarines 
and  underwater  objects,  and  low  frequency  systems  detect  large  objects  over  further 
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distances  [24],  An  example  of  SONAR  imagery  is  Fig.  4c  which  shows  the  caldera  of  an 
underwater  volcano  in  the  Lau  Basin  near  the  Samoan  Islands  [23].  The  bottom  portion  of 
Fig.  4c  is  imagery  available  via  the  Google  Maps™  mapping  service  and  is  a  combination 
of  SONAR  data  and  satellite  altimetry,  a  RADAR  technique.  SONAR  can  also  be  an 
example  of  a  passive  RS  system  when  used  in  “listen”  mode  only;  detecting  naturally 
occurring  sounds  rather  than  the  reflection  of  artificially  generated  sound  waves  [20], 

More  common  passive  RS  techniques  include  thermal,  panchromatic,  multispectral, 
and  polarimetric  imaging.  Thermal  detectors  such  as  thermocouples  and  bolometers  sense 
emitted  IR  radiation.  Any  minute  change  in  temperature  triggers  a  change  to  the  resistivity 
in  the  detector’s  circuitry  in  the  case  of  a  bolometer  and  a  change  in  the  voltage  across  two 
dissimilar  metal  leads  in  the  case  of  a  thermocouple  [51].  These  changes  in  resistivity  or 
voltage  can  be  processed  into  an  image  such  as  that  shown  in  Fig.  5a  displaying  energy 
loss  through  a  window  header  [41].  Thermal  imaging  applications  include  detection  of 
energy  radiated  from  water  vapor  at  the  top  of  the  atmosphere  revealing  weather  patterns, 
military  technology  such  as  night  vision  goggles,  and  satellites  that  can  detect  the  launch 
of  ballistic  missiles.  The  main  advantage  of  thermal  imaging  is  its  day  and  night  capability, 
but  one  shortcoming  is  the  slow  response  time  of  most  thermal  detectors  compared  to 
photodetectors. 


(a)  Thermal  [41]  (b)  Multispectral  [9]  (c)  Polarimetric  [47] 

Figure  5:  Passive  Remote  Sensing  Product  Examples 
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Surrounding  the  thermal  image  in  the  center  of  Fig.  5a  is  a  panchromatic  image. 
Panchromatic  (or  grayscale)  images  are  formed  by  electronically  converting  the  number  of 
photons  incident  onto  a  detector  over  a  broad  wavelength  spectrum,  usually  in  the  visible 
band,  into  a  value  representing  the  radiance  of  the  scene  of  interest.  This  type  of  image  only 
conveys  information  about  the  apparent  brightness  of  the  scene  and  is  unable  to  convey 
color  information.  Visible  multispectral  imagery  is  the  process  of  taking  panchromatic 
images  over  multiple  narrow  wavelength  bands  of  interest.  Each  spectral  image  appears 
as  grayscale  until  processing  techniques  are  used  to  combine  the  images  into  a  “true”  or 
“false”  color  image.  The  most  common  mutispectral  bands  are  red,  green,  and  blue  which 
are  combined  to  produce  a  true  color  image  showing  a  scene  as  it  might  look  to  the  human 
eye.  Figure  5b  shows  a  false  color  multispectral  image  of  the  1988  Yellowstone  fires  taken 
by  a  twelve  band  imager  aboard  a  National  Aeronautics  and  Space  Administration  (NASA) 
owned  aircraft  helping  authorities  track  the  fire’s  progression  through  the  park  [6, 8,9]. 

2.1.2  Passive  Multispectral  Polarimetric  Remote  Sensing. 

The  least  intuitive  RS  technique  is  undoubtedly  polarimetry.  Polarization  is  a  property 
of  all  TEM  waves,  so  its  applications  are  not  limited  to  any  specific  wavelength  region. 
Most  RS  textbooks  only  discuss  polarization  with  regards  to  its  application  in  RADAR 
systems,  including  a  text  solely  committed  to  the  topic  of  applications  of  polarization  in 
remote  sensing  [5,21,26].  A  recent  text  by  Schott  [36]  recognizes  that  “lack  of  data 
continues  to  prevent  many  application  scientists  from  evaluating  the  value  of  PI  remote 
sensing. . .  it  is  the  need  for  better  understanding  by  the  potential  user  community  of  the 
subtleties  of  PI  that  prompted  the  development  and  structuring  of  this  text.”  Passive  Pi’s 
utility  lies  in  the  fact  that  when  natural  (unpolarized)  light  scatters  off  a  natural  or  man¬ 
made  object  its  polarization  state  is  altered.  Figure  5c  is  a  polarimetric  image  of  an  outdoor 
scene  containing  a  thick  line  of  trees,  a  grassy  area,  and  two  pickup  trucks  parked  in  the 
shade  under  the  trees.  The  interested  reader  is  referred  to  Tyo’s  [47]  Fig.  1  for  true  color  and 
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long-wave  IR  intensity  representations  of  the  scene.  As  shown  in  Tyo’s  work,  the  trucks 
are  not  visible  in  the  true  color  or  the  long- wave  IR  intensity  images.  The  polarimetric 
image,  however,  distinctly  reveals  the  vehicles  in  the  upper  middle  portion  of  the  image. 
As  shown  in  this  work,  PI  can  unearth  more  than  simply  the  presence  of  man-made  targets. 
If  MS  functionality  is  added  to  a  polarimetric  imager,  this  work’s  proposed  algorithm  could 
recover  the  target’s  refractive  index,  thus  determining  its  surface  material.  Even  in  the 
absence  of  polarization  content  (which  is  typical  at  near-normal  incidence)  the  true  color 
images  produced  by  a  MS  polarimetric  imager  would  still  be  useful  for  target  identification. 

2.2  Title  Components  Theory 

With  the  immense  diversity  of  possible  RS  techniques  established  and  the  case  made 
for  combining  multispectral  and  polarimetric  RS,  the  focus  of  Chapter  II  is  narrowed  to 
detailing  the  theory  behind  each  component  of  this  work’s  title. 

2.2.1  Material  Characterization. 

As  introduced  in  section  1.2,  the  permittivity,  conductance,  and  permeability  of  a 
material  are  conveyed  by  its  index  of  refraction.  Characterizing  a  material  involves 
finding  a  physical  theory  linking  the  refractive  index  to  measurable  parameters  of  interest 
allowing  for  either  direct  inversion  to  recover  the  refractive  index  or  solving  via  least 
squares  methods.  In  this  work  the  index  of  refraction  is  found  using  Matrix  Laboratory 
(MATLAB)’s  Isqcurvefit  function,  a  non-linear  data-fitting  least  squares  solving  method. 
In  order  for  Isqcurvefit  to  recover  a  material’s  refractive  index,  the  material  must  first  be 
classified.  Classifying  a  material  as  either  dielectric  or  metallic  is  less  complicated  than  full 
characterization  in  that  precisely  solving  for  the  refractive  index  is  not  necessary.  Given 
an  empirical  data  set,  a  researcher  can  decide  on  some  threshold  criteria  with  which  to 
classify  the  material.  Full  characterization  via  RS  is  desirable  for  intelligence  gathering, 
non-destructive  inspection,  and  many  other  situations  when  physical  contact  with  a  material 
is  denied  or  unsuitable. 
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2.2.2  Multispectral. 

Having  defined  multispectral  imagery  in  section  2.1.1,  this  section  focuses  on  one 
complication  applicable  to  this  work  that  arises  when  taking  multiple  narrow  band  images 
of  a  scene:  chromatic  aberration.  Figure  6  shows  dispersion  (discussed  further  in 
section  2.3.2)  of  incident  light  rays  causing  the  focal  point  of  each  wavelength  to  differ  [30], 
The  optical  components  of  this  work’s  experimental  setup  were  installed  to  optimize 
performance  at  the  633nm  focal  point,  thus  imaging  any  other  wavelength  produces  a 
defocused  image.  Without  a  complete  system  redesign,  this  problem  could  have  been 
remedied  by  finding  the  system’s  Point  Spread  Function  (PSF)  at  each  wavelength  and 
deconvolving  each  defocused  image,  but  a  simpler  approach  (manual  focusing)  was  used 
in  the  interest  of  time. 


Figure  6:  Chromatic  Aberration  Caused  by  a  Convex  Lens  [30] 

2.2.3  Polarimetric. 

Figure  7  depicts  how  EM  (light)  waves  propagate;  the  electric  field  is  perpendicular  to 
the  magnetic  field  and  both  are  perpendicular  to  the  direction  of  propagation  (£,  77,  and  z 
in  the  figure,  respectively).  The  polarimetric  state  of  the  wave  shown  in  Fig.  7  is  considered 
linear  and  vertical  since  the  projection  of  the  electric  field  vector  as  it  propagates  traces  out 
a  vertical  line  in  the  xy  plane.  The  wave  in  Fig.  8a  is  also  considered  linear,  but  is  oriented 
at  45°  with  respect  to  the  x  axis  as  shown  by  the  dotted  trace  of  its  <5V  and  &y  components 
onto  the  xy  plane  [36],  Figure  8  also  shows  how  a  TEM  wave’s  polarization  state  may  be 
circular  or  elliptical  and  can  appear  to  rotate  counterclockwise  as  in  Fig.  8b  or  clockwise 
as  in  Fig.  8c  [36].  Although  not  perfectly  so,  sunlight  is  considered  unpolarized  meaning 
the  electric  field  is  randomly  oriented  as  the  wave  propagates. 
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Figure  7:  Fluctuations  in  the  Electric  and  Magnetic  Fields  Give  Rise  to  TEM  Wave 
Propagation  in  the  z  Direction  [36] 


Figure  8:  Polarization  States  [36] 

Polarization  state  cannot  be  measured  by  any  sensor  directly,  but  can  be  calculated 
from  the  electric  field’s  irradiance.  Introduced  by  G.  G.  Stokes  in  1852,  Stokes  parameters 
are  used  to  represent  light’s  polarization  state  [16].  An  interpretation  of  how  the  parameters 
characterize  polarization  state  is  offered  in  Table  1  [36],  Conventionally,  each  Stokes 
parameter  is  normalized  to  the  S0  value,  allowing  better  representation  of  how  much  each 
parameter  contributes  to  the  wave’s  polarimetric  state. 


Table  1:  Meaning  of  Stokes  Parameter  Values,  Modified  from  Schott  [36] 


Parameter  Positive 

So  Represents  total  S  field  irradiance 

S\  Horizontal  polarization  dominates 

S2  +45°  polarization  dominates 

St,  Right-hand  circular  dominates 


Negative 

N/A 

Vertical  polarization  dominates 
-45°  polarization  dominates 
Left-hand  circular  dominates 
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The  last  Stokes  parameter,  S3,  is  not  applicable  to  this  work.  When  illuminated  by 
a  natural  source,  light  scattered  from  an  object  has  been  shown  to  contain  an  insignificant 
amount  of  circularly  polarized  light  [7, 18,22,37,49].  As  a  consequence  of  this  observation, 
the  polarization  analyzer  used  by  this  work  need  only  consist  of  a  linear  polarizer.  The 
mathematical  implication  of  this  simplification  is  explained  in  section  2.3.3. 

2.2.4  Imagery. 

In  the  context  of  this  research,  images  are  a  series  of  grayscale  irradiance 
measurements  collected  of  an  unknown  target  using  a  multispectral  polarimeter.  The 
system  required  calibration  and  characterization  before  any  experimental  data  was 
collected.  If  the  system  was  ideal,  its  resolution  would  be  limited  only  by  diffraction 
caused  by  its  entrance  pupil,  a  circle.  The  diffraction  pattern  measured  when  imaging  a 
point  source  is  the  system’s  PSF.  In  the  frequency  domain,  the  inverse  Fourier  transform 
of  the  PSF  is  the  system’s  Optical  Transfer  Function  (OTF)  and  conveys  which  frequencies 
are  attenuated  or  passed  by  the  system.  The  Modulation  Transfer  Function  (MTF)  is  the 
absolute  value  of  the  OTF  and  is  typically  normalized  and  plotted  as  a  function  of  cutoff 
frequency.  One  way  to  measure  the  quality  of  an  imaging  system,  thus  aiding  in  calibration, 
is  to  determine  its  MTF.  Any  operational  optical  system  will  suffer  from  aberrations 
(defocus,  spherical,  etc.)  that  degrade  image  quality.  The  aberrations  are  apparent  in  the 
system’s  MTF.  A  set  of  orthogonal  polynomials  called  Zernikes  can  be  used  to  determine 
how  many  waves  of  each  aberration  exist  in  a  system.  Each  aberration’s  polynomial  is 
added  together  and  applied  as  a  phase  to  a  complex  pupil  function.  The  pupil  function  is 
then  Fourier  transformed  to  obtain  the  amplitude  spread  function.  The  absolute  value  of 
the  amplitude  spread  function  gives  the  PSF,  and  then  as  before,  the  MTF  is  found.  The 
mathematical  form  of  diffraction  limited  MTF  curves  are  known  and  are  plotted  along  with 
a  system’s  MTF  to  show  how  near  an  imaging  system  is  to  diffraction  limited  performance. 
Data  showing  the  performance  of  the  instrument  used  in  this  work  is  given  in  section  3.1.2. 
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2.3  Mathematical  Framework 


With  the  theoretical  context  of  this  work  defined,  this  section  presents  the  specific 
mathematical  framework  used  to  produce  this  work’s  results. 

2. 3. 1  Complex  Index  of  Refraction. 

A  material’s  complex  index  of  refraction  q  =  n  -  ]k  is  a  combination  of  two  optical 
constants:  index  of  refraction  n  and  extinction  coefficient  k.  The  refractive  index  n  is 
defined  as  the  ratio  of  phase  velocity  of  light  in  vacuum  to  the  phase  velocity  of  light  in  the 
medium  ( n  =  c/v).  The  extinction  coefficient  k  is  related  to  the  exponential  decay  of  the 
wave  as  it  passes  through  the  medium  [2]  and  can  be  thought  of  as  representing  absorption. 

2.3.2  Dispersion  Equations. 

Since  the  underlying  physical  parameters  that  define  the  refractive  index  (conductivity, 
permeability,  and  permittivity)  depend  on  wavelength,  so  does  q  itself.  Dispersion  refers 
to  the  phenomena  whereby  a  medium’s  index  of  refraction  is  frequency  (wavelength) 
dependent  [16].  Evidence  of  dispersion  is  shown  by  shining  a  light  containing  all  visible 
wavelengths  (a  white  light)  onto  a  prism  as  seen  in  Fig.  9  [10].  A  material’s  refractive  index 
is  measured  over  many  wavelengths  to  produce  its  dispersion  curve  such  as  that  for  gold  in 
Fig.  10  [28], 


16 


Figure  10:  Empirical  Dispersion  Curve  of  Gold,  h(-),  k(-  -)  [28] 


Empirical  dispersion  curves  are  available  for  a  wide  variety  of  materials  and  are 
approximated  by  dispersion  equations  such  as  the  Cauchy  equation: 

A; 


the  Sellmeier  equation: 


and  the  Drude  model: 
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where  a>  =  2nc/A  ,  cop  is  plasma  frequency,  and  r  is  electron  relaxation  rate  [13, 16]. 
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The  constants  A0,  A,,  and  C,  in  Eq.  (1)  and  Eq.  (2)  are  borne  out  of  what’s  happening 
to  a  dielectric  at  the  atomic  level  as  it’s  bombarded  by  a  harmonic  EM  wave.  The  Lorentz 
model  of  light/matter  interaction  likens  the  bonds  between  electrons  and  their  nuclei  to  the 
classical  mechanics  concept  of  a  small  mass  connected  by  a  spring  to  a  larger  mass.  An 
incident  electric  field  induces  displacement  of  the  electrons  which  then  undergo  a  spring¬ 
like  restoring  force  toward  equilibrium.  The  natural  resonance  frequencies  in  a  volume 
of  molecules  due  to  these  spring-like  oscillations  give  rise  to  multiple  oscillatory  modes. 
The  dispersion  equation  constants  give  the  wavelength  at  which  those  modes  occur  and 
the  relative  strength  of  each  mode.  The  interested  reader  is  referred  to  Hecht  [16]  section 
3.5.1  for  more  detail.  Although  the  Lorentz  model  was  not  published  until  the  late  1890s, 
Cauchy  published  his  equation  for  dispersion  in  1835  based  on  his  own  theory  (or  arguably 
Fresnel’s  theory)  of  how  particles  responded  to  an  external  force  [3],  Sellmeier’s  equation 
better  predicts  dispersion  over  a  wider  wavelength  range  and  was  published  in  1871  [2], 
Since  this  work  is  limited  to  the  visible  portion  of  the  EM  spectrum,  the  Cauchy  equation 
is  sufficient  for  estimating  the  refractive  index  of  dielectrics. 

As  mentioned  above,  the  Lorentz  model  treats  electrons  as  bound  to  the  nucleus  by  a 
harmonic  restoring  force  thus  representing  a  non-conductive  dielectric.  The  Drude  model 
is  identical  to  the  Lorentz  model,  but  sets  the  restoring  force  to  zero  thus  allowing  the 
electrons  to  move  freely  representing  a  metal  [2].  Free  electrons  and  positive  ions  within 
a  metal  may  be  thought  of  as  a  plasma  whose  density  oscillates  at  a  natural  frequency  a>p, 
the  plasma  frequency,  which  for  most  metals  occurs  within  the  visible  or  Ultraviolet  (UV) 
bands  [16],  Asa  consequence  of  the  Drude  model,  the  refractive  index  of  metals  is  found 
using  Eq.  (3).  This  work’s  application  of  dispersion  equations  is  explained  in  section  3.2.1. 

Every  material  has  wavelength  bands  within  their  dispersion  curve  that  follow  what 
are  called  “normal”  and  “anomalous”  dispersion.  Normal  dispersion  describes  a  band  over 
which  the  real  part  of  the  refractive  index  decreases  with  increasing  wavelength,  whereas 
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anomalous  dispersion  occurs  when  n  increases  with  wavelength  [16].  Where  these  regions 
occur  within  the  EM  spectrum  depend  on  where  strong  absorption  occurs  for  that  material. 
Absorption  occurs  when  the  incident  light  has  a  high  enough  energy  to  excite  electrons 
within  the  material  from  a  lower  energy  state  to  a  higher  energy  state,  thus  “absorbing” 
the  incoming  energy.  Most  metals  experience  no  absorption  within  the  visible  band  and 
instead  absorb  in  the  UV  region,  meaning  they  undergo  anomalous  dispersion  in  the  visible 
and  reflect  almost  all  visible  light.  As  shown  in  Fig.  11  [17],  the  combination  of  the  three 
additive  colors  red,  green,  and  blue  give  most  metals  a  white  or  silver  appearance.  Gold  and 
copper  are  examples  of  metals  that  experience  absorption  in  the  visible  band.  Gold  absorbs 
incident  blue  light,  thus  explaining  how  the  remaining  visible  light  that  either  reflects  off 
of  or  is  re-emitted  by  gold  (red  plus  green)  gives  the  material  its  yellow  appearance.  This 
also  means  gold  experiences  normal  dispersion  in  the  visible  band.  As  a  consequence 
of  this  information  and  the  fact  that  this  work’s  data  is  collected  in  the  visible  band, 
it  is  expected  that  upon  evaluation  of  the  metallic  samples  used  in  this  work  (gold  and 
aluminum),  gold’s  refractive  index  values  will  decrease  with  wavelength  while  aluminum’s 
n  values  will  increase  with  wavelength. 


Figure  11:  RGB  Model  of  Additive  Light  [17] 
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2.3.3  Stokes  Vectors  and  Mueller  Calculus. 


The  meaning  of  each  Stokes  parameter  was  given  in  section  2.2.3.  Using  the  axes 
convention  given  in  Fig.  8  it  is  clear  at  any  given  time  t  the  electric  field  vector  6  has 
amplitude  components  Sqx  and  &oy  in  the  x  and  y  planes,  respectively.  Together  these 
components  form  a  transverse  harmonic  plane  wave  such  as  that  expressed  by 

<5  =  Sx(z,  t)  +  Sy(z,  t ) 

Sx(z,  t)  =  <5q  v  cos (kz  -  tot)  i  (4) 

Sy(z,  t )  =  S0y  cos  (kz  -  a>t  +  e)  j 

where  z  is  the  position  of  the  wave,  e  is  the  relative  phase  difference  between  the  wave’s 
components  (for  linear  polarization  e  =  0,  n,  or  2n),  k  is  the  wavenumber  (not  to  be 
confused  with  the  extinction  coefficient),  and  a>  is  the  wave’s  angular  frequency.  The  most 
general  state  of  polarization  is  elliptical  with  linear  and  circular  polarization  being  special 
cases  thereof.  In  a  derivation  shown  by  Hecht  [16],  the  time  and  position  dependencies  can 
be  removed  from  Eq.  (4)  to  express  the  electric  field  vector’s  projected  ellipse  on  the  xy 
plane  given  by 


After  some  algebraic  manipulation  including  time  averaging  as  Schott  mentions  is  shown 
in  Stokes  [36, 39],  the  projected  ellipse  can  be  expressed  in  the  form 

(Si  +  &i)2  =  (Si  -  Si)2  +  (2&1&1  cos  6)2  +  (2&ox&oy  Sin  e)2  (6) 


The  Stokes  parameters  for  a  plane  wave  are  thus  defined  in  Schott  [36]  as 

50  =  &l  +  &l 

51  =  &l  -  So, 

Si  =  2  &1&1  cos  6 
S3  —  28{lx8ih:  sin  € 


(7) 
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The  parameters  depend  only  on  measurable  irradiances  and  are  conventionally  combined 
into  a  Stokes  vector  as 

S  =  [S0  Sl  S2]t  (8) 


The  Stokes  vector  given  in  Eq.  (8)  contains  only  So  through  S2  because,  as  explained  in 
section  2.2.3,  this  work  is  only  concerned  with  linearly  polarized  light  and  S3  represents  the 
amount  of  circularly  polarized  light.  Examples  of  Stokes  vectors  representing  commonly 
produced  linear  polarimetric  states  include 


[1  0 

or 

Unpolarized 

[1 

i  or 

Horizontal 

[1  -1 

or 

Vertical 

[1 

0  1]T 

Linear  45° 

which  agree  with  the  information  presented  in  Table  1 .  The  Degree  of  Linear  Polarization 
(DOLP)  is  the  degree  to  which  light  is  linearly  polarized  (a  value  from  0  representing  no 
linear  polarization  to  1  meaning  light  is  fully  linearly  polarized)  and  is  expressed  in  terms 
of  the  Stokes  parameters  as 


DOLP 


+  S‘ 


(10) 


The  Angle  of  Polarization  (AOP)  is  the  angle  the  electric  field  vector  makes  with  respect  to 
the  x  axis  and  is  defined  in  terms  of  Stokes  parameters  as 


AOP  =  ^  arctan  j  (1 1) 

The  goal  of  passive  PI  is  to  calculate  an  “input”  Stokes  vector  S)n  which  gives  the 
polarization  state  of  light  scattered  off  a  target  of  interest.  The  standard  equation  used  to 
calculate  Sm  is 

Sout  =  Msys  •  Sin  (12) 


where  the  S0  parameter  of  the  “output”  Stokes  vector  Sout  is  the  total  irradiance  measured 
by  the  PI  system  (Sjjx  +  Sjjy)  and  Msys  is  the  PI  system’s  Mueller  matrix.  Mueller  matrices 
characterize  the  diattenuation  and  retardance  of  polarization  elements  and  are  known  for 
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many  optical  components  [2],  Any  component  of  a  PI  system  that  may  affect  Sjn  before  it 
reaches  the  imager  such  as  a  polarizer  or  retarder  must  be  accounted  for  in  the  system 
Mueller  matrix.  If  a  system  has  more  than  one  polarization  altering  component,  each 
component’s  individual  Mueller  matrix  can  be  multiplied  to  form  one  overall  system 
Mueller  matrix.  Further  discussed  in  section  3.1.2,  calculation  of  5in  is  straightforward 
with  multiple  measurements  of  S  [jut  and  a  known  Msys. 

2.3.4  Fresnel  Reflectance. 

The  equations  governing  how  light  interacts  with  matter  are  a  good  place  to  start 
when  searching  for  a  physical  model  relating  the  measurable  Stokes  parameter  S  °ut  to 
the  refractive  index.  Using  Maxwell’s  equations,  the  law  of  reflection,  and  Snell’s  law 
of  refraction,  it  can  be  shown  [18]  that  the  Fresnel  reflection  coefficients  given  by 

_Ucos6i-A)2  +  B2' 

■'  “  [(cos 0,  +A)2  +  52 

(A  -  sin  6i  tan  d,)2  +  B2 
p  (A  +  sin  6i  tan  O,)2  +  B2  s 


C  =  4n2k2  +  D2 
D  =  n2  -  k2  -  sin2  d; 


are  dependent  on  a  material’s  refractive  index.  To  understand  the  connection  between  the 
Fresnel  reflection  coefficients  Rs  and  Rp  and  the  Stokes  parameters,  it  helps  to  visualize 
the  light/matter  interaction.  Figure  12  shows  Hyde’s  [18]  representation  of  two  cases  of 
polarized  light  incident  on  a  smooth  medium  of  refractive  index  rj  =  n-jk  where  E,  H,  and 
k  correspond  to  the  wave’s  electric  field,  magnetic  field,  and  propagation  direction  vectors 
respectively,  and  the  superscripts  i,  r,  and  t  correspond  the  to  the  incident,  reflected,  and 
transmitted  portions  of  the  wave.  The  specular  plane  in  Fig.  12  (and  in  general)  is  defined 
as  the  plane  containing  the  material’s  surface  normal  and  the  incident  and  reflected  light’s 
propagation  direction  vectors.  If  the  entire  electric  field  vector  lies  perpendicular  to  the 
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(a)  S-Pol  (b)  P-Pol 

Figure  12:  S-Pol  and  P-Pol  Scattering  Geometries  [18] 


specular  plane  as  in  Fig.  12a,  the  incident  light  is  referred  to  as  perpendicular  polarization 
or  S-Pol.  If  the  entire  electric  field  vector  is  in  the  specular  plane  as  in  Fig.  12b,  the  incident 
light  is  referred  to  as  parallel  polarization  or  P-Pol  [18].  Derivations  in  Torrance  [46]  and 
He  [15]  give  the  sought  after  physical  model  relating  the  Fresnel  reflectance  equations  to 
the  Stokes  parameters  as  applied  to  the  case  of  specular  reflection  from  a  smooth  surface: 


DOP  = 


Rs  ~  Rp 

Rs  +  Rp 


(14) 


2.3.5  Bidirectional  Reflectance  Distribution  Function. 

The  previous  section  applied  to  the  case  of  reflection  off  a  smooth  surface.  A 
surface  looks  smooth  to  incident  light  if  the  surface  roughness  is  much  smaller  than  the 
light’s  wavelength.  In  the  visible  band,  400  -  700nm,  most  surfaces  look  “rough”  so 
the  model  stated  above  must  be  reconsidered.  Figure  13  introduces  the  concept  of  a 
Bidirectional  Reflectance  Distribution  Function  (BRDF).  Instead  of  modeling  reflection 
as  purely  specular  as  with  the  image  on  the  far  left,  a  BRDF  models  the  reflection  of  light 
more  realistically;  somewhere  between  specular  and  diffuse  reflection.  A  BRDF  describes, 
mathematically,  how  reflected  light  is  directionally  distributed  within  the  hemisphere  above 
the  material’s  surface.  When  constructing  a  BRDF,  whether  empirically  or  analytically, 
a  material  surface  geometry  (smooth/rough/otherwise)  must  be  known  or  assumed.  The 
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Near-perfect 
diffuse  reflector 


diffuse  reflector 
("Lambertian  surface") 


Figure  13:  Light  Scattering  Models  [21] 


most  commonly  used  model  is  to  think  of  a  surface  as  a  composition  of  microfacets  (similar 
to  the  facets  of  a  diamond,  but  on  a  much  smaller  scale),  each  an  ideal  specular  reflector 
obeying  Snell’s  law  (0,  =  6r)  with  their  own  Fresnel  reflectivity  and  their  own  geometry  with 
respect  to  the  light  source  and  RS  receiver  [29].  The  BRDF  surface  and  viewing  geometry 
used  in  this  work  mimics  Hyde’s  [18]  and  is  given  in  Fig.  14. 


24 


The  directional  distribution  of  a  TEM  wave’s  polarimetric  content  can  equivalently  be 
described  by  a  Polarimetric  Bidirectional  Reflectance  Distribution  Function  (pBRDF).  The 


pBRDF  used  in  this  work  is  given  by 


fjiiflu  0r,  4>) 


l  exp  [- tan2  0/(2cr2)] 

To  \7a  2\  4~n  A  A  m/°i,  Or,  <p) 

(2^)(4cr2)  cos4  0  cos  0r  cos  8j 


(15) 


where  /;7  is  the  element  in  the  /h  row,  /lh  column  of  the  Mueller  matrix  pBRDF,  m/7  is  the 
element  in  the  jth  row,  Ith  column  of  the  Mueller  matrix  for  reflection,  </>  is  the  azimuthal 
angle,  cr  is  the  surface  roughness  parameter,  and  6  is  the  microfacet  angle  of  orientation 
relative  to  the  mean  surface  normal  given  by 


6  =  cos  1 


cos  +  cos  6r 
2  cos/? 


(16) 


where 


cos  2/3  =  cos  0r  cos  +  sin  6r  sin  dt  cos  </> 


(17) 


[18,29,42],  Since  this  work  is  confined  to  the  specular  plane  where  <p  =  n,6  and /?  simplify 
to 

9  =  (0,.  -  d,)/2  /3  =  (0r  +  9,)/2  (18) 


As  shown  in  Hyde  [18]  and  Thilak  [42],  using  the  pBRDF  equation  to  represent  DOP 
once  again  yields  Eq.  (14).  Substituting  in  Rs  and  Rp  from  (13)  and  simplifying  produces 
the  equation  underlying  all  of  this  work’s  results: 


P(n,k,/3 )  = 


2 A  sin2/?  cos/? 


A2  cos 2  {3  +  sin4/?  +  B 2  cos2/? 


VC  +  D 


B 


C  =  4  n2k2  +  D 


VC-D 


(19) 


D  =  rr  -  k2  -  sin2/? 


where  P  is  the  DOP,  or  in  the  case  of  this  work,  DOLP. 
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III.  Methodology 


With  the  relevant  theory  and  mathematical  framework  in  place,  this  chapter  describes 
the  method  used  to  apply  those  concepts  to  the  problem  statement,  thus  producing 
a  hypothesis.  Experimentation  (rather  than  simulation)  was  used  to  test  this  work’s 
hypothesis,  so  a  description  of  the  experimental  instrument,  analysis  of  the  instrument’s 
optical  and  polarimetric  performance,  and  how  the  instrument  was  used  to  collect  data  is 
discussed.  The  goal  of  this  work  is  not  only  to  classify,  but  to  characterize  unknown  targets, 
so  the  methods  used  for  classification  and  characterization  are  then  discussed. 

3.1  Experimental  Implementation 

Data  was  collected  using  the  Multispectral  Stokes  Imager  (MS  SI)  instrument  (Fig.  15) 
in  the  Optical  Turbulence  Estimation,  Compensation,  and  Simulation  (OPTECS)  laboratory 
at  the  Air  Force  Institute  of  Technology.  A  description  of  the  instrument,  its  performance, 
and  how  it  was  used  to  collect  data  follows. 

3.1.1  Instrument  Description. 

The  light  source  (Fig.  16a)  is  a  mounted,  collimated,  and  diffused  warm  white  Light 
Emitting  Diode  (LED)  attached  to  a  rotating  base  which  allows  for  modification  of  the 
angles  of  incidence  and  reflection.  An  optical  filter  mount  holds  the  target  in  place  in 
a  rotation  stage.  Actuation  of  the  rotating  source  base  and  the  target  rotation  stage  is 
controlled  via  a  software  interface.  The  light  scattered  off  the  target  is  reflected  off  a 
mirror,  which  sends  it  down  the  rest  of  the  optical  train.  The  Atmospheric  Turbulence 
Simulator  (ATS)  portion  of  the  instrument  (Fig.  16b)  was  not  used  for  this  research,  but 
was  left  in  place  for  future  extensions  to  this  and  other  OPTECS  experiments.  Following 
the  ATS,  the  beam  is  split  with  one  half  the  light  proceeding  to  a  point  polarimeter  and 
the  other  half  proceeding  through  to  the  Stokes  imager  (Fig.  16c).  The  point  polarimeter 
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automatically  calculates  the  Stokes  parameters  and  DOLP  (among  many  other  variables) 
allowing  for  verification  of  input  Stokes  “truth”  states  and  measured  DOLP  values.  The 
MSSI  consists  of  a  rotating  linear  polarizer,  a  manual  filter  wheel  with  six  lOnm  bandgap 
spectral  filters  ranging  in  wavelength  from  450  -  700nm  in  50nm  increments,  a  focusing 
lens,  an  aperture  to  test  diffraction  limited  performance  and  reduce  irradiance  as  necessary, 
and  finally,  the  camera.  Given  the  camera’s  focal  plane  size  of  2560  x  1920  pixels,  a  pixel 
pitch  of  2.2pm,  and  a  Field  of  View  (FOV)  of  2.125in.  (0.053975m),  the  MSSI  demagnifies 
the  target  to  10.38  percent  of  its  original  size  meaning  each  pixel  images  roughly  21pm  of 
the  target. 


Figure  15:  Experimental  Instrument 
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(a)  Source  and  Target 


(b)  ATS 


(c)  Analyzer  and  Imaging  Optics 

Figure  16:  Experimental  Instrument  Components 
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3.1.2  Instrument  Performance. 

First,  the  system’s  imaging  quality  as  compared  to  the  diffraction  limited  case  is 
considered,  followed  by  an  analysis  of  the  polarimetric  response  of  the  system.  An  imaging 
system’s  performance  is  often  expressed  by  its  MTF  as  compared  to  an  ideal,  diffraction 
limited  system.  Computation  of  a  system’s  MTF  was  discussed  in  2.2.4.  To  measure  the 
MSSI’s  MTF,  a  15pm  pinhole  was  inserted  into  the  target  holder  and  illuminated  with 
a  collimated  laser  beam.  The  image  of  this  approximated  point  source  is  defined  as  the 
system’s  PSF  and  MATLAB  was  used  to  perform  the  necessary  mathematical  steps  to  arrive 
at  the  MTF.  The  MSSI’s  measured  MTF  is  shown  as  the  dashed  blue  curve  in  Fig.  17. 


Figure  17:  System  Defocused,  Focused  Corrected,  and  Theoretical  MTF  Curves 
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The  diffraction  limited  MTF  is  shown  as  the  green  line  in  Fig.  17  and  was  calculated 


using 


MTF(f)  = 


2 

cos  1 

(X) 

\  _  (X)2 

71 

(2/0) 

2/0  V 

A  \2f0)  J 

0 

if  /  <  2/o 

otherwise 


(20) 


where  /  =  +  fy,  fx  =  x/Azu  fy  =  y/Azi,  and  /0  is  the  system’s  cutoff  frequency 

defined  as  /o  =  d/2Azi  [14].  The  MTF  was  found  using  the  laser’s  wavelength  A  =  635nm, 
a  circular  pupil  diameter  of  size  d  =  12.3mm,  and  an  image  distance  of  Zi  =  10.5cm. 

The  MTF  is  a  function  of  spatial  frequency  /  and  the  closer  the  system’s  MTF  curve 
is  to  the  diffraction  limited  case,  the  higher  the  system’s  spatial  frequency  resolution 
will  be,  signifying  a  finer  resolution  imaging  capability.  Since  the  MSSI’s  MTF  curve 
was  not  diffraction  limited,  an  attempt  was  made  to  quantify  the  aberrations  causing  the 
instrument’s  sub-optimal  performance.  Potential  aberrations  were  modeled  using  Zemike 
polynomials.  The  plot  of  the  Zemike  estimated  MTF  and  the  approximated  aberration 
values  is  shown  as  the  red  line  in  Fig.  17.  Clearly  the  Zernike  estimated  MTF  does  not 
match  the  system  MTF  perfectly,  but  any  change  to  the  given  polynomial  coefficient  values 
caused  an  increase  in  mean  squared  error.  Therefore,  the  given  estimation  was  deemed 
sufficient.  The  MSSI’s  only  correctable  system  aberration  is  defocus.  The  change  in 
optical  path  length  necessary  to  account  for  the  system’s  defocus  was  calculated  to  be 
-0.14068pm.  Measuring  to  that  precision  is  impractical  for  this  application,  so  an  Air 
Force  resolution  bar  chart  was  inserted  into  the  target  holder  and  the  system  was  visually 
focused.  The  shape  of  the  camera’s  pixel  intensity  histogram  distinctly  fluctuates  as  the 
image  reaches  and  then  loses  focus.  When  the  shape  of  the  Air  Force  bar  chart’s  histogram 
indicated  the  best  possible  focus,  the  pinhole  was  re-inserted  as  the  target  and  another  PSF 
measurement  was  taken  revealing  marked  improvement  as  shown  by  the  black  curve  in 
Fig.  17.  The  tendency  of  the  MTF  to  “overshoot”  the  diffraction  limit  is  termed  apodization 
and  is  caused  by  any  introduction  of  absorption  into  the  pupil  [14]. 
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With  the  instrument’s  imaging  quality  analyzed,  its  polarization  response  was  sought. 
As  mentioned  in  section  2.3.3,  a  Mueller  matrix  accounts  for  the  polarization  altering 
effects  of  components  within  the  optical  system.  The  MSSI’s  only  polarization  component 
is  a  rotating  linear  polarizer,  which  has  a  known  Mueller  matrix,  given  by 


Ma  = 


1 

2 


1 

cos  26 
sin  26 


cos  26 
cos  22 6 
sin  26  cos  26 


sin  26 

sin  26  cos  26 
sin  22  6 


(21) 


where  6  is  the  angle  of  the  linear  polarizer  (further  referred  to  as  the  analyzer)  with  respect 
to  the  x  axis.  Notice  a  Mueller  matrix  is  nominally  a  4  x  4  matrix,  but  since  this  work  does 
not  deal  with  circular  polarization  the  fourth  row  and  fourth  column  are  unnecessary.  If  the 
analyzer  was  the  only  element  contributing  to  the  system’s  polarization  response,  then  the 
input  Stokes  vector  could  be  found  by  simply  solving 

Sout  =  Ma  •  Sin  (22) 

given  a  measured  output  Stokes  vector.  Unfortunately,  passive  PI  systems  are  only  capable 
of  directly  measuring  the  total  irradiance  of  the  field  incident  upon  the  polarimeter  (given 
by  the  first  parameter  S  [j111  of  the  output  Stokes  vector).  The  other  parameters  of  the 
output  Stokes  vector,  S°ut  and  .S’™'1,  are  not  measurable  making  direct  inversion  of  Eq.  (22) 
impossible.  Due  to  this  complication  another  method  is  sought  for  finding  Sin. 

A  set  of  “truth”  input  Stokes  vectors  was  generated  by  a  precision  Meadowlark  linear 
polarizer  inserted  directly  after  a  633nm  LED  source  (Fig.  18).  Using  the  known  input 
Stokes  vectors,  a  measured  set  of  Sgut  parameters,  and  the  data  reduction  method  described 
below,  it  was  determined  that  the  system’s  polarization  response  was  not  caused  by  the 
analyzer  alone.  Knowing  there  must  be  some  other  component(s)  affecting  the  polarization 
state,  Eq.  (22)  was  updated  to  include  an  unknown  system  Mueller  matrix  Msys: 

Sout  =  Ma  •  Msys  •  Sin  (23) 
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Figure  18:  Meadowlark  Linear  Polarizer  and  633nm  LED  source 


Recovering  the  unknown  system  Mueller  matrix  and  future  unknown  input  Stokes 
vectors  is  possible  using  the  data  reduction  technique  described  by  Bass  [2]  and  Tyo  [47], 
The  only  vector  and  matrix  elements  needed  from  Eq.  (23)  to  express  S  JJut  are 
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where  Msys  is  broken  down  into  generic  Mueller  matrix  elements.  To  demonstrate  the 
simple  linear  nature  of  this  method,  Eq.  (24)  is  multiplied  out  as 

SoUt(0)  =  1(6)  =^[5o(lm0o  +  cos  2  6  ■  mw  +  sin  2  6  ■  m2  o) 

+  ^“(lmoi  +  cos  2 6  ■  mn  +  sin  26  •  m2i)  (25) 

+  ^“(lmfn  +  cos  20  •  m\2  +  sin  26  •  m22)] 


A  system  of  linear  equations  can  be  built  up  by  collecting  multiple  irradiance  S  )jut 
measurements  at  different  analyzer  angles  8j  as  in 
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The  S  QUt(0/)  total  irradiance  values  of  Eq.  (26)  are  measured  by  the  camera,  the  6j 
values  of  the  analyzer  Mueller  matrix  are  known,  and  the  input  Stokes  vector  values  are 
known  given  the  truth  state  generated  by  the  Meadowlark  linear  polarizer.  Knowing  these 
values  it  is  possible  to  determine  all  nine  Msys  elements  with  a  minimum  of  nine  linear 
equations.  An  over-determined  system  approach  was  taken  in  this  work.  In  total,  16  linear 
equations  were  produced  by  collecting  16  irradiance  data  sets  in  the  form  of  videos  which 
were  later  separated  into  S  image  frames  and  averaged  over  each  set  of  frames.  The  16 
videos  correspond  to  the  following  algorithm  of  data  collection: 

1.  Set  up  and  turn  on  the  warm  white  diode  source 

2.  Set  Meadowlark  linear  polarizer  (input  Stokes  vector)  to  horizontal  state 

3.  Rotate  analyzer  to  6j  =  0°,  45°,  90°,  and  135°  taking  a  video  at  each  angle 

4.  Repeat  steps  2  and  3  for  remaining  input  Stokes  states  (vertical,  +45°,  and  -45°) 

The  analyzer  angles  used  are  arbitrary,  but  the  angles  chosen  for  this  work  are 

conventionally  used  and  correspond  to  the  angles  used  in  the  fabrication  of  Kulkarni  and 
Gruev’s  [19]  new  division  of  focal  plane  multispectral  polarimeter,  a  sensor  the  MSSI  could 
be  adapted  to  use  in  the  future.  MATLAB  was  used  to  set  up  the  matrices  in  Eq.  (26)  and 
use  of  the  proper  matrix  operations  yielded  the  unknown  Msys.  Following  this  procedure 
the  calculated  Msys  was  inserted  back  into  Eq.  (23)  and  the  correct  “truth”  input  Stokes 
vectors  were  recovered,  proving  the  validity  of  the  calculated  Msys.  Since  polarization 
response  is  wavelength  dependent,  this  procedure  was  repeated  for  each  of  the  six  spectral 
filters  used  in  this  work.  The  only  complication  encountered  during  data  acquisition  at 
each  wavelength  was  chromatic  aberration,  as  discussed  in  section  2.2.2.  The  presence 
of  chromatic  aberration  required  the  system  be  manually  focused  after  each  spectral  filter 
change.  With  Msys  known  for  each  wavelength,  the  only  remaining  experimental  task  was 
to  collect  polarimetric  data  of  targets  at  multiple  wavelengths. 
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3.1.3  Multispectral  Polarimetric  Data  Collection  Algorithm. 

The  algorithm  used  to  collect  target  data  closely  resembles  that  used  to  find  Msys  for 
each  wavelength  and  is  as  follows. 

1.  Choose  target  and  insert  into  the  target  holder 

2.  Rotate  target  and  source  to  achieve  desired  specular  angle  of  incidence 

3.  Turn  on  source 

4.  Rotate  analyzer  to  90°  to  allow  maximum  irradiance  through  system 

5.  Rotate  filter  wheel  to  spectral  filter  with  highest  transmittance 

6.  Adjust  source  intensity  to  ensure  no  camera  pixels  are  saturated 

7.  Rotate  filter  wheel  to  desired  measurement  wavelength  filter 

8.  Manually  focus  system 

9.  Center  target  area  of  interest  using  software  cross  hairs  and  target  surface  feature 

10.  Rotate  analyzer  to  6  =  0°,  45°,  90°,  and  135°  taking  a  set  of  frames  at  each  angle 

11.  Repeat  steps  7  through  10  for  each  measurement  wavelength 

12.  Repeat  steps  2  through  1 1  for  each  desired  angle  of  incidence 

This  manual,  iterative  process  of  data  collection  is  time  consuming  and  fraught 
with  misregistration  issues  (discussed  further  in  subsequent  sections),  but  suffices  for  the 
purpose  of  demonstrating  this  work’s  novel  data  analysis  algorithm.  Potential  steps  to 
automate  and  improve  upon  the  data  collection  algorithm  are  discussed  in  section  5.2. 

3.2  Data  Analysis  Algorithm 

Since  the  variable  of  interest  to  this  work,  the  index  of  refraction,  is  different  at  each 
wavelength  (as  evidenced  by  dispersion),  the  data  analysis  algorithm  must  consist  of  more 
than  arranging  the  collected  data  into  a  set  of  linear  DOLP  equations  and  attempting  to 
solve  for  //.  In  order  to  create  a  well-posed  set  of  equations,  the  non-linear  least  squares 
solver  used  to  find  rj  (MATLAB’s  Isqcurvefit  function)  must  be  solving  for  a  variable  that 
does  not  vary  with  wavelength.  Given  the  theoretical  background  and  the  mathematical 


34 


framework  from  Chapter  II,  a  novel  hypothesis  was  formed  that  the  dispersion  equations 
for  dielectric  and  metallic  materials  could  be  substituted  in  place  of  n  and  k  in  the  pB RDF- 
based  DOLP  equation  (Eq.  (19)).  The  non-linear  least  squares  solver  could  then  search 
for  the  wavelength-independent  dispersion  equation  constants  (A,  for  dielectrics  and  ojp  for 
metals)  that  satisfied  the  empirical  DOLP  values  on  the  left-hand  side  of  Eq.  (19).  In  order 
to  test  this  theory,  the  material  must  first  be  classified  as  dielectric  or  metallic. 

3.2.1  Material  Classification. 

In  his  work,  Hyde  [18]  uses  a  pBRDF  model  to  simulate  the  DOLP  for  different 
materials  of  differing  surface  roughness.  He  notes  that  the  maximum  DOLP  values 
achieved  for  dielectric  and  metallic  materials  are  approximately  1  and  0.3,  respectively. 
He  also  mentions  that  other  material  classification  techniques  [4, 45,  52]  take  advantage 
of  the  generalization  that  dielectric  materials  strongly  polarize  scattered  light,  whereas 
metallic  materials  weakly  polarize  scattered  light.  This  work  uses  the  maximum  DOLP 
values  found  by  Hyde  as  simple  threshold  classifiers.  In  order  to  use  the  threshold  values 
to  classify  each  pixel,  the  DOLP  of  each  pixel  must  first  be  found.  The  frames  of  irradiance 
data  of  a  target  from  the  data  collection  algorithm  are  averaged  and  used  to  find  the  input 
Stokes  vector  of  each  pixel  at  each  wavelength.  The  DOLP  of  each  pixel  is  found  for  each 
wavelength  using  the  calculated  Stokes  parameters  as  in  Eq.  (10).  Then,  the  DOLP  values 
of  one  pixel  across  all  wavelengths  are  averaged.  This  one,  wavelength  averaged  DOLP 
value  per  pixel  determines  whether  that  pixel  is  classified  as  a  dielectric  or  a  metal.  If  the 
averaged  DOLP  is  >  0.3,  the  pixel  is  classified  as  dielectric  and  if  DOLP  is  <  0.3,  the  pixel 
is  classified  as  metal.  Hyde  warns,  though,  that  DOLP  varies  over  angle  of  incidence,  so 
the  model  should  not  be  construed  as  a  rule,  but  rather  a  guide.  At  near-normal  incident 
angles  a  better  classification  technique  may  be  to  use  the  derivative  of  the  DOLP  values 
over  wavelength  due  to  the  tendency  of  all  materials  to  polarize  scattered  light  very  weakly 
in  that  regime.  The  derivative  would  be  near  zero  for  metals  and  measurable  for  dielectrics. 
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3.2.2  Material  Characterization. 

This  work’s  data  analysis  code  has  two  non-linear  least  squares  solving  routines,  one 
for  dielectrics  that  substitutes  the  Cauchy  dispersion  equation  into  the  n  variable  of  Eq.  (19), 
and  one  for  metals  that  substitutes  the  Drude  dispersion  equation  into  n  and  k  of  Eq.  (19). 
Once  a  pixel  is  classified,  the  measured  DOLP  values  of  that  pixel  for  each  wavelength  are 
sent  to  the  appropriate  least  squares  solver  routine  and  the  Cauchy  coefficients  or  the  plasma 
frequency  is  found.  The  solver  returns  the  Cauchy  coefficients  or  the  plasma  frequency 
which  are  then  plugged  back  into  the  appropriate  dispersion  equation  to  recover  n  and 
k.  The  entire  classification  and  characterization  process  is  repeated  for  each  pixel  and  is 
depicted  by  Fig.  19.  The  non-linear  least  squares  solving  routines  both  use  MATLAB’s 
Isqcurvefit  function  to  find  the  dispersion  equation  constants.  Reasonable  bounds  were 
placed  on  the  potential  constant  values  of  between  0  and  5  for  the  Cauchy  coefficients  [28] 
and  1  x  1014  and  1  x  1017  for  the  plasma  frequency  [13]. 
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Figure  19:  Proposed  Per-Pixel  MS  PI  Data  Analysis  Algorithm  Flow  Chart 
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IV.  Results  and  Discussion 


As  mentioned  in  section  2.1.2,  the  benefits  of  combining  MS  and  PI  include  their 
individual  conventional  abilities  to  produce  true  color  images  and  to  distinguish  between 
man-made  and  natural  objects,  respectively.  A  set  of  data  products  is  presented  below  to 
demonstrate  these  benefits.  As  shown  by  this  work’s  novel  algorithm,  another  benefit  of 
combining  MS  and  PI  is  the  ability  to  classify  and  characterize  an  unknown  target .  Figures 
showing  classification  and  characterization  results  and  tables  displaying  measured  DOLP 
and  refractive  index  values  are  given  below  to  show  the  utility  of  this  work’s  algorithm. 

4.1  Multispectral  and  Polarimetric  Imaging  Data  Products 

Much  like  the  Mueller  matrices  that  accounted  for  the  system’s  polarization  response 
at  each  wavelength,  a  set  of  multiplication  factors  was  necessary  to  account  for  the  spectral 
response  of  the  filters  (Fig.  20b  shows  one  example)  and  the  camera  (Fig.  20c).  The 
multiplication  factors  normalize  the  contribution  of  each  wavelength  to  ensure  the  true 
color  image  accurately  represents  the  target  and  were  calculated  as  the  inverse  of  each 
filter’s  transmission  percentage  multiplied  by  the  camera’s  quantum  efficiency.  Accounting 
for  the  spectral  response  of  the  system  is  not  necessary  when  calculating  the  Stokes 
parameters  since  those  are  based  on  irradiance  differences  measured  at  the  same  wavelength 
rather  than  across  multiple  wavelengths.  Furthermore,  the  DOLP  values  are  calculated 
using  normalized  Stokes  parameters,  so  any  compensation  for  spectral  response  would 
vanish  during  the  normalization  process. 

Each  wavelength’s  S  [jut  image  was  multiplied  by  the  appropriate  factor  and  then 
normalized  to  the  highest  irradiance  value  across  all  wavelengths.  The  normalized  values 
were  then  assigned  to  a  Red,  Green,  and  Blue  (RGB)  image  data  matrix  and  with 
MATLAB’s  image  function  the  true  color  images  were  produced  (Fig.  21).  MATLAB’s 
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image  function  takes  an  n  x  m  x  3  matrix  of  RGB  values  and  specifies  the  true  color  of  each 
pixel  from  a  color  map.  The  abundance  of  650nm  light  incident  on  the  camera,  particularly 
when  combined  with  the  high  reflectivity  of  a  metallic  target,  limited  the  source’s  output 
intensity  and  caused  a  general  lack  of  dielectric  brightness  in  the  mixed  material  targets. 


Wavelength  (nm) 


(a)  Source  Spectrum  [44] 


Wavelength  (nm) 


(b)  650nm  Spectral  Filter  Transmission  [43] 


(c)  Camera  Quantum  Efficiency  [1] 

Figure  20:  OPTECS  MSSI  System  Spectral  Characteristics 
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(e)  Aluminum  and  Green  Paint 


(f)  Circuit  Board 


Figure  21:  Experimental  Multispectral  True  Color  Images 


39 


As  expected,  glass  (Fig.  21a)  shows  as  transparent  (interpreted  by  MATLAB  as  white, 
a  combination  of  all  colors).  Some  pixel  misregistration  can  be  seen  where  the  magenta, 
cyan,  and  yellow  spots  do  not  line  up  to  form  white.  A  set  of  natural  and  man-made  leaves 
(Fig.  21b)  are  mainly  green  with  some  blue  and  pink  displaying  misregistration.  A  sample 
of  LabSphere  InfraGold  has  hints  of  its  characteristic  yellow  coloration  (Fig.  21c),  but 
misregistration  is  so  abundant  in  this  image  that  interpretation  of  the  target  as  gold  solely 
from  the  true  color  image  is  difficult.  Data  for  the  rough  Aluminum  (Al)  target  shown  in 
Fig.  21d  was  collected  at  only  450nm,  550nm,  and  650nm.  Unfortunately,  the  approximate 
wavelengths  necessary  to  construct  a  true  color  image  align  better  with  450nm,  550nm, 
and  700nm  (the  filters  used  for  all  other  data  collection),  so  the  rough  Al  target  has  a 
misleading  orange/pink  hue.  A  target  composed  of  Al  with  applied  green  paint  clearly 
shows  the  Al  region  and  a  less  distinct  green  tint  to  the  paint  region  due  to  the  limitation  on 
source  intensity  mentioned  earlier.  The  gold  fingers  and  epoxy  solder  mask  of  the  circuit 
board  (Fig.  2 If)  are  identifiable  although  the  color  of  the  solder  mask  is  not  consistent  with 
the  original  target.  The  apparent  blue  color  of  the  green  solder  mask  could  be  due  to  the 
imprecise  nature  of  manually  focusing  the  instrument  or  lack  of  camera  dynamic  range. 

4.2  Classification  Results 

Figure  22  shows  the  classification  result  for  each  target,  with  DOLP  >  0.3  =  dielectric 
and  DOLP  <  0.3  =  metal.  The  only  surprising  result  was  that  of  the  set  of  leaves  (Fig.  22b). 
The  classification  is  correct,  but  the  DOLP  values  of  the  natural  leaf  are  much  higher  than 
those  of  the  man-made  leaf,  which  is  the  opposite  of  the  expected  result.  The  issue  here 
lies  not  in  the  theory  regarding  DOLP  values  for  man-made  versus  natural  targets,  but  in 
the  data  collection  scenario.  For  an  operational  remote  sensing  system  the  target  would  be 
much  further  away  from  the  sensor  than  is  possible  to  duplicate  in  a  lab  environment.  As 
a  result,  when  an  operational  system  images  a  scene  it  would  “see”  a  canopy  of  tree  leaves 
which  diffusely  scatters  fight  rather  than  seeing  only  one  highly  specular  reflecting  leaf. 
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(e)  Aluminum  and  Green  Paint 


(f)  Circuit  Board 


Figure  22:  Target  DOLP  Averaged  over  All  Wavelengths  (left)  and  Material  Classification 
Results  (right) 


Classifying  individual  pixels  of  a  mixed  material  target  can  lead  to  misclassification 
as  shown  by  the  circuit  board  and  the  A1  and  green  paint  targets  in  Fig.  22.  Averaging  the 
DOLP  over  an  adjacent  grouping  of  pixels  and  then  classifying  could  increase  the  accuracy 
of  the  classification.  For  instance,  the  small  area  of  misclassified  material  of  the  A1  and 
green  paint  target  would  be  classified  correctly  if  averaged  with  the  surrounding  (lower) 
DOLP  values.  It  may  be  undesirable  to  average  over  many  pixels  depending  on  the  system’s 
resolution,  but  in  this  case,  with  one  pixel  imaging  only  2 1  pm  of  the  target,  the  loss  of  a  bit 
of  resolution  for  more  accurate  classification  would  be  an  acceptable  trade-off. 
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Section  1.2  mentions  the  importance  of  correct  pixel  registration  between  irradiance 
images  over  all  wavelengths  to  this  work.  One  drawback  to  the  manual  nature  of  the  MSSI 
is  pixel  misregistration.  Throughout  the  data  collection  process  an  obvious  feature  on 
the  target  surface  was  chosen  and  the  cross  hairs  of  the  camera  software  were  used  to 
re-center  the  target’s  area  of  interest  (typically  200  x  200  pixels)  on  that  feature.  Any 
misalignment  suffered  between  wavelength  measurements  is  visible  in  the  normalized  S  o 
images  of  each  target  (two  examples  are  shown  in  Fig.  23).  The  misalignment  in  the  case 
of  the  leaves  target  is  easy  to  see  when  examining  the  center  (550  nm)  image  as  opposed 
to  the  outer  images.  The  consequence  of  this  misalignment  is  reduction  in  classification 
and  characterization  accuracy.  When  the  DOLP  values  of  a  pixel  from  each  wavelength 
are  passed  to  the  classification  algorithm  they  are  averaged.  If  one  value  is  inflated  due 
to  misregistration  it  could  cause  misclassification  of  the  pixel.  Similarly,  with  respect  to 
the  non-linear  least  squares  solver,  the  DOLP  values  used  to  determine  the  constants  in  the 
dispersion  equations  would  not  represent  an  optimal  curve  fit  scenario  if  one  or  more  of 
them  was  abnormal,  thus  potentially  causing  mischaracterization  of  that  pixel. 


(b)  Leaves 

Figure  23:  Normalized  So  Parameter  Showing  700nm  (left),  550nm  (center),  and  450nm 
(right)  Content  of  Gold  and  Leaves 
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4.3  Characterization  Results 


As  mentioned  in  section  3.2.2,  once  a  pixel  is  classified  its  set  of  DOLP  values 
(one  for  each  measurement  wavelength)  is  sent  to  the  Cauchy  or  Drude  model  function 
as  appropriate.  The  metallic  plasma  frequency  or  set  of  dielectric  Cauchy  coefficients  is 
found  via  non-linear  least  squares  fit  to  the  DOLP  values.  Then  the  coefficients  or  plasma 
frequency  are  substituted  into  the  appropriate  dispersion  equation  to  recover  n  and  k.  Below 
are  the  characterization  results  and  analysis  for  each  target  examined  by  this  work. 

Calculation  of  a  theoretical  value  for  DOLP  at  each  wavelength  is  possible  using 
Eq.  (19)  along  with  the  published  t]  values  for  each  target  material  cited  throughout  this 
section  and  the  known  angle  of  incidence.  Table  2  gives  the  results  of  those  calculations 
as  well  as  the  experimentally  obtained  DOLP  values  averaged  over  wavelength.  The 
experimental  DOLP  values  of  targets  composed  of  only  one  material  were  calculated 
directly  by  averaging  over  every  pixel  of  every  wavelength.  Mixed  material  targets  were 
either  visually  examined  to  determine  an  approximate  average  DOLP  value  or  segmented 
using  a  threshold  DOLP  value  and  then  averaged  over  each  segment.  As  expected, 
the  dielectric  materials  exhibit  a  much  higher  DOLP  than  the  metallic  materials.  Also 
noteworthy  is  the  accuracy  of  the  experimental  values  despite  only  defining  the  “macro” 
incident  angle  rather  than  the  more  precise  angle  of  incidence  that  occurs  on  the  micro¬ 
facet  level.  DOLP  accuracy  is  key  to  the  least-squares  solver’s  ability  to  correctly  identify 
the  Cauchy  coefficients  and  plasma  frequency  per  pixel  ultimately  leading  to  the  correct  rj. 

Table  2:  Theoretical  and  Experimental  DOLP  Values  of  each  Target,  Averaged  over  A 


Theoretical  Experimental 


Material 

6i  =  6r  =  45° 

di  =  dr= 60° 

9i  =  0r  =  45° 

ei=er=60° 

Glass 

0.82 

0.98 

0.83 

0.94 

Natural  Leaf 

0.83 

0.98 

0.77 

0.91 

Gold 

0.17 

0.21 

0.14 

0.14 

Aluminum 

0.03 

0.06 

0.07 

0.10 

Green  Paint 

0.80 

0.94 

0.83 

0.95 

Epoxy 

0.83 

0.98 

0.80 

0.90 
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4.3.1  Dielectric  Characterization. 

The  Cauchy  equation,  as  discussed  in  section  2.3.2,  is  a  good  approximation  of 
dielectric  dispersion  in  the  visible  band.  A  limitation  of  this  model,  however,  is  its 
applicability  to  only  non-lossy  dielectrics  (non-conductive  materials  that  are  transparent  in 
the  visible  band).  Depending  on  the  accuracy  desired  and  the  number  of  MS  bands  available 
on  a  given  RS  platform,  the  Sellmeier  equation  could  be  used  along  with  a  minimum  of 
2 i  +  1  (index  value  i,  not  imaginary  i)  wavelength  measurements  to  recover  n  and  k.  This 
work  used  the  Cauchy  equation  for  dielectric  pixels  and  as  a  result  only  recovered  n  while 
setting  k  =  0.  Figure  24  shows  the  experimental  n  values  for  the  glass  target  and  the  target 
consisting  of  a  set  of  leaves.  The  mean  experimental  n  value  at  each  wavelength  for  the 
glass  target  is  in  close  agreement  with  the  known  n  values  of  BK7  glass  (Table  3)  computed 
using  the  Sellmeier  coefficients  listed  on  Schott’s  N-BK7  datasheet  [35], 


(b)  Leaves 

Figure  24:  Dielectric  Target  n  Values  at  450nm  (left),  550nm  (center),  and  700nm  (right) 
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Table  3:  Mean,  Standard  Deviation,  and  Percent  Error  of  n  for  Glass 


/t(nm) 

E[n] 

O'  n 

Known  n 

%  Error 

450 

1.4267 

0.0202 

1.52532 

6.47 

550 

1.4260 

0.0206 

1.51852 

6.09 

700 

1.4255 

0.0208 

1.51306 

5.79 

A  set  of  refractive  index  values  is  available  for  natural  foliage  [48]  and  is  given  in 
Table  4.  Judging  by  the  n  values  indicated  by  the  color  scale  in  Fig.  24,  the  natural  leaf 
does  appear  to  have  a  refractive  index  near  1.5  over  each  wavelength. 


Table  4:  Refractive  Index  of  Natural  Foliage  [48] 


/t(nm) 

Known  n 

450 

1.527 

550 

1.509 

700 

1.493 

4.3.2  Metal  Characterization. 

The  usefulness  of  the  Drude  model  depends  on  the  ability  to  break  rj  from  Eq.  (3)  into 
an  expression  for  the  real  part  n  and  the  extinction  coefficient  k  and  then  substituting  those 
expressions  back  into  Eq.  (19).  The  refractive  index  is  defined  as  the  square  root  of  the 
complex  relative  dielectric  constant  (er  =  s\  -)s2)  where  n  and  k  are  related  to  s\  and  s2  in 


Palik  [28]  by 


n 


(s2+s2)  +  si 


and  k 


l(e?  +  £l)  ~  £i 


(27) 


and  £\  and  s2  are  defined  in  Palik  [28]  as 


S\  =  n  -  k2  =  1 


2  2 
cozpT2 


and  s2  =  2  nk  = 


oj2pt 


(28) 


(1  +  oj2t2)  a>(l  +  co2t2) 

The  electron  relaxation  rate  r  (representing  the  time  between  collisions)  is  dependent  on 
frequency  and  can  be  described  by 


ojz 


T  = 


4n(TOJ2 


(29) 
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where  cr  is  the  material’s  dc  conductivity  [28].  In  general,  r  increases  with  wavelength  and 
is  cited  to  equal  approximately  10“ 14  for  most  metals  [2, 13, 16].  Use  of  the  correct  value 
for  r  is  very  important  to  accurately  calculating  rj.  A  good  deal  of  effort  was  put  forth 
to  find  consistent  published  values  of  r  for  gold  and  aluminum,  but  with  values  ranging 
from  10“ 15  to  10-12  (often  without  a  wavelength  reference)  another  method  of  finding  r 
was  sought.  By  manually  forcing  the  least-squares  solver  to  use  gold’s  plasma  frequency 
top  of  13.6973  x  1015  rad/s  [13],  it  was  possible  to  find  values  of  r  that  produced  acceptably 
accurate  n  and  k  values  for  gold  (r  =  3.25  x  10-14,  1.25  x  10-15,  and  4  x  10-15  for  450,  550, 
and  700nm  respectively).  The  functional  form  of  r\  used  by  the  least  squares  solver  should 
ideally  not  contain  any  wavelength-dependent  variables.  For  this  reason,  every  instance  of  r 
in  Eq.  (28)  should  be  replaced  by  its  own  functional  form  (Eq.  (29))  and  the  dc  conductivity 
cr  should  be  treated  as  another  constant  for  the  least  squares  solver  to  find. 

Table  5  and  Fig.  25  show  the  characterization  results  of  the  gold  target  with  Q,  =  0,  = 
45°  and  known  n  and  k  values  from  Palik  [28],  The  results  shown  for  gold  use  the  r  values 
generated  by  forcing  the  plasma  frequency  to  equal  its  known  value.  When  the  raw  data 
was  evaluated  for  presentation  the  bounds  on  ojp  were  relaxed  leaving  the  solver  to  work  as 
intended.  Allowing  the  solver  to  determine  ojp  (rather  than  forcing  the  known  a>p)  produced 
a  mean  plasma  frequency  of  9.62  x  1015  rad/s  representing  an  error  of  29.7%.  Since  gold 
absorbs  near  450nm,  the  n  value  at  that  wavelength  is  especially  sensitive  to  errors  in  r, 
possibly  explaining  the  relatively  larger  percent  error  of  n  at  450nm.  Also,  it  is  possible  the 
n  value  at  450nm  reported  in  Palik  is  not  quite  precise  itself. 


Table  5:  Mean,  Standard  Deviation,  and  Percent  Error  of  n  and  k  for  the  Gold  Target 


/t(nm) 

EM 

&~n 

Known  n 

%  Error,, 

E[k] 

o-k 

Known  k 

%  Error* 

450 

0.7423 

0.1926 

1.5023 

50.59 

1.7044 

0.7016 

1.8785 

9.27 

550 

0.3421 

0.1005 

0.3593 

4.79 

2.5539 

0.9554 

2.6913 

5.11 

700 

0.1725 

0.0536 

0.1611 

7.08 

3.4094 

1.2160 

3.9521 

13.73 
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(a)  Gold  n 


Figure  25:  Experimental  Refractive  Index  Component  Values  at  450nm  (left),  550nm  (center), 
and  700nm  (right)  for  the  Gold  Target 


The  r  values  used  for  the  gold  data  were  used  again  when  evaluating  the  rough  A1  data. 
Figure  26  shows  the  resulting  n  and  k  images  for  each  wavelength.  The  mean  experimental 
n  and  k  values  given  in  Table  6  are  compared  to  those  given  by  Palik  [28].  The  data 
furnished  in  Palik  were  acquired  for  bulk  polished  Al,  so  the  errors  may  be  a  product  of 
differing  surface  roughness  and  the  inaccuracy  of  r.  The  mean  plasma  frequency  calculated 
by  the  least  squares  solver  was  2.33  x  10“ 16  rad/s,  representing  only  a  3%  error  from  the 
known  value  of  2.4  x  10“16  rad/s  [13].  Given  a  nearly  ideal  <op,  the  refractive  index  values 
likely  would  also  have  been  nearly  ideal  had  the  correct  (but  unknown)  r  values  been  used. 


Table  6:  Mean,  Standard  Deviation,  and  Percent  Error  of  n  and  k  for  the  Rough  Al  Target 

/t(nm) 

EM 

Known  n 

%  Error,, 

E[k] 

o-k 

Known  k 

%  Error* 

450 

1.5412 

0.9567 

0.6332 

132.4 

4.4067 

3.0504 

5.4544 

19.2 

550 

0.7489 

0.4775 

1.0159 

26.28 

6.2781 

4.2457 

6.6273 

5.27 

650 

0.3345 

0.2153 

1.5580 

78.53 

7.5875 

5.0864 

7.7124 

1.62 
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(b)  Rough  A1  k 

Figure  26:  Experimental  Refractive  Index  Component  Values  at  450nm  (left),  550nm  (center), 
and  650nm  (right)  for  the  Rough  A1  Target 


4.3.3  Mixed  Material  Characterization. 

The  final  test  of  this  work’s  material  characterization  algorithm  was  to  characterize 
targets  composed  of  dielectric  and  metallic  materials.  Shown  and  discussed  below  are  the 
results  for  the  A1  with  green  paint  and  the  circuit  board  targets. 

The  Drude  model  values  of  r  were  originally  found  for  a  normally  dispersive  metal 
(gold),  so  it  is  no  surprise  that  the  experimental  n  values  for  the  A1  region  in  Fig.  27  decrease 
rather  than  increase  as  they  should.  The  decreasing  values  of  n  may  lead  an  observer  to 
think  the  target  was  copper  or  gold,  but  one  look  at  the  true  color  image  (if  collected  at  the 
correct  wavelengths)  would  reveal  the  metal  to  be  colorless,  invalidating  the  n  values.  If  the 
least  squares  solver  had  been  designed  to  find  r  the  n  values  would  be  more  accurate.  The 
green  paint  region  closely  matches  the  n  values  in  Table  7,  but  due  to  this  work’s  use  of  the 
Cauchy  rather  than  the  Sellmeier  equation,  the  k  values  of  green  paint  were  not  recovered. 
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(b)  A1  with  Green  Paint  k 

Figure  27:  Experimental  Refractive  Index  Component  Values  at  450nm  (left),  550nm  (center), 
and  700nm  (right)  for  the  A1  with  Green  Paint  Target 


Table  7:  Known  Refractive  Index  of  Green  Paint,  Extrapolated  from  Data  in  Shell  [37] 


T(nm) 

Extrapolated  n 

Extrapolated  k 

450 

1.46285 

0.296025 

550 

1.45815 

0.298675 

700 

1.4511 

0.30265 

Taking  a  look  at  the  circuit  board  k  value  images  in  Fig.  28b  (in  which  the  limits 
of  the  data  were  manipulated  only  to  show  better  contrast),  it  is  easy  to  spot  the  regions 
that  were  likely  misclassified.  The  boundaries  of  the  gold  fingers  are  clearly  defined  and 
are  filled  with  a  high  density  of  large  k  values.  The  areas  of  diffuse  pixels  were  likely 
misclassified.  The  misclassification,  and  corresponding  mischaracterization,  could  be  a 
result  of  the  semi-transparent  nature  of  the  epoxy  in  between  the  fingers.  The  cause  could 
be  misregistration  if  at  one  wavelength  the  gold  fingers  overlapped  an  area  of  imaged  epoxy 
at  another  wavelength,  but  judging  from  the  alignment  shown  by  the  normalized  So  images 
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in  Fig.  29,  misregistration  is  likely  not  the  cause.  Another  cause  of  the  misclassification 
could  be  scattering  occurring  from  a  certain  depth  within  the  circuit  board,  essentially 
displaying  the  depth  of  the  gold  finger,  but  this  is  unlikely.  Large  portions  of  the  k  values  of 
the  gold  at  each  wavelength  agree  with  the  known  values  for  gold  given  in  Table  5.  Also, 
the  n  values  of  the  gold  finger  regions  decrease  with  wavelength,  suggesting  they  represent 
gold  or  copper.  The  true  color  image  in  Fig.  21f  confirms  a  gold-colored  appearance. 

The  paths  of  the  copper  traces  leading  from  the  gold  regions  are  fairly  visible  and  have 
a  higher  refractive  index  value  than  their  immediate  surroundings.  Since  the  real  part  of 
copper’s  refractive  index  is  a  bit  lower  than  gold’s  in  the  visible  band  (1.15-0.2  [28]),  there 
may  be  a  highly  insular  material  protecting  the  copper  traces  and  the  areas  in  between  the 
gold  fingers,  producing  higher  n  values  in  those  regions.  The  n  values  elsewhere  are  typical 
of  epoxy  used  as  a  solder  mask  on  printed  circuit  boards  [11]. 


(a)  Circuit  Board  n 


-0.2  -0.15  -0.1  -0.05  0  0.05  0.1  0.15  0.2  -0.2  -0.15  -0.1  -0.05  0  0.05  0.1  0.15  0.2  -0.2  -0.15  -0.1  -0.05  0  0.05  0.1  0.15  0.2 


(b)  Circuit  Board  k 

Figure  28:  Experimental  Refractive  Index  Component  Values  at  450nm  (left),  550nm  (center), 
and  700nm  (right)  for  the  Circuit  Board  Target 
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Figure  29:  Normalized  So  Parameter  Showing  Red  (left),  Green  (center),  and  Blue  (right) 
Content  of  Circuit  Board 


Throughout  this  chapter,  surface  irregularities  have  been  apparent  in  each  target’s  n 
and  k  images  (edge  of  where  the  green  paint  meets  the  aluminum,  the  fabric-like  woven 
texture  of  the  man-made  leaf,  roughness  of  the  aluminum,  etc.).  Surface  features  such  as 
these  scatter  more  light  parallel  to  the  surface  feature  than  perpendicular  to  the  feature, 
increasing  the  DOLP  in  those  regions  and  affecting  the  n  and  k  values. 
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V.  Conclusions  and  Future  Work 


5.1  Conclusions 

This  thesis  proved  the  author’s  hypothesis  correct  within  acceptable  error  limits  and 
in  a  lab  environment.  With  known  6t  and  6,  angles  and  DOLP  measurements  collected  at 
multiple  wavelengths,  the  target’s  index  of  refraction  can  be  estimated  via  a  non-linear  least 
squares  fit  using  two  dispersion  equations.  Each  pixel  of  an  imaged  target  was  classified  as 
dielectric  or  metallic  and  was  subsequently  characterized  using  a  non-linear  least  squares 
solver  and  .  A  pBRDF-derived  equation  for  the  DOLP  was  given  each  pixel’s  set  of 
DOLP  values  (one  for  each  measurement  wavelength)  and  the  least  squares  solver  found 
the  dispersion  equation  constants  that  generated  the  best  fit  to  the  experimental  data.  The 
dispersion  equation  constants  were  substituted  back  into  the  dispersion  equations  yielding 
the  complex  index  of  refraction  for  each  pixel. 

Although  the  characterization  results  are  impressive,  the  most  promising  results  are 
displayed  in  Table  2  comparing  the  theoretical  and  experimental  DOLP  values.  As  long 
as  the  experimental  DOLP  values  are  reasonably  close  to  the  theoretical  values  the  least 
squares  solver  has  a  good  chance  of  finding  dispersion  equation  constants  that  produce 
accurate  refractive  index  values.  The  modifications  to  and  extensions  of  this  work  described 
in  the  next  section  will  only  serve  to  increase  the  algorithm’s  accuracy. 

The  combination  of  multispectral  and  polarimetric  imagery  is  powerful.  The  true 
color  images  correlated  with  the  refractive  index  data  certainly  has  the  potential  to  increase 
target  identification  accuracy.  Even  in  the  absence  of  polarization  content,  as  is  the  case 
for  near-normal  collection  geometries,  the  true  color  image  alone  would  provide  useful 
target  identification  data.  While  many  new  remote  sensing  data  interpretation  techniques 
add  layers  of  ever-increasing  complexity,  this  novel  algorithm  provides  elegant  simplicity 
and  proven  potential  to  the  task  of  high  value  target  detection  and  identification. 
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5.2  Future  Work 


5.2.1  Material  Classification  and  Characterization  Algorithm  Improvements. 

Rather  than  using  a  simple  thresholding  method  to  classify  each  pixel,  the  technique 

developed  by  Wolff  [52]  may  lead  to  better  classification  results.  Using  the  Sellmeier 
in  place  of  the  Cauchy  dispersion  equation  to  characterize  dielectrics  would  improve  the 
accuracy  of  the  characterization  algorithm  as  it  would  allow  for  the  recovery  of  k.  Also,  as 
mentioned  multiple  times,  replacing  r  in  the  Drude  model  dispersion  equation  and  solving 
for  the  wavelength-independent  dc  conductivity  would  ensure  more  accurate  r  values  and 
likely  much  more  accurate  refractive  index  values. 

5.2.2  Expanding  Algorithm  Region  of  Validity. 

The  results  presented  in  this  work  display  the  material  characterization  algorithm’s 
validity  at  an  incident  angle  of  6t  =  45°.  Collection  of  additional  data  at  many  other  angles 
of  incidence  may  prove  or  disprove  the  usefulness  of  this  algorithm  over  a  wide  range  of 
collection  angles.  The  previously  stated  limitation  of  little  polarization  content  at  near¬ 
normal  geometries  could  be  quantified  by  collecting  data  as  close  to  the  surface  normal 
vector  as  the  instrument  allows  (approximately  6,  =  17°). 

5.2.3  Accounting  for  Microfacet  Geometry. 

The  AOP  can  be  used  to  obtain  the  microfacet  normal  for  each  pixel  or  a  certain  sized 
group  of  pixels.  To  the  author’s  knowledge  a  mature  and  concise  technique  to  accomplish 
this  does  not  yet  exist,  but  if  accomplished  could  provide  a  more  accurate  value  of  [3  for  the 
DOLP  equation  and  would  allow  expansion  of  this  algorithm  to  out-of-plane  scattering.  In 
this  work  (3  reduces  to  simply  6,  due  to  the  condition  of  specular  reflection,  but  this  would 
not  be  the  case  for  out-of-plane  scattering.  The  interested  reader  can  find  more  information 
regarding  AOP  techniques  and  the  full / 3  equation  in  Hyde  [18],  Robles  [33],  and  Stolz  [40]. 
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5.2.4  Incorporation  of  Atmospheric  Turbulence. 

Many  RS  systems  such  as  astronomic  telescopes  on  Earth  and  imaging  systems  used 
from  aircraft  are  affected  by  atmospheric  turbulence.  Temperature  and  pressure  differences 
in  the  atmosphere  cause  pockets  of  air  with  differing  indices  of  refraction,  distorting  a  RS 
system’s  acquired  image.  The  algorithm  developed  here  could  be  incorporated  into  the 
work  accomplished  by  Hyde  [18]  to  create  an  algorithm  capable  of  producing  turbulence- 
corrected  MS  polarimetric  images. 

5.2.5  System  Automation. 

Future  instantiations  of  the  OPTECS  Lab  MSSI  instrument  would  greatly  benefit  from 
automation.  The  possible  misregistration  issues  mentioned  in  section  4.2  could  be  resolved 
by  the  use  of  a  calibrated  automatic  filter  wheel  and  an  auto-focusing  camera.  Automation 
would  also  save  data  collection  time,  allowing  more  samples  to  be  tested.  Replacing 
the  filter  wheel,  analyzer,  and  the  camera  with  a  division  of  focal  plane  MS  polarimeter 
such  as  that  produced  by  Kulkarni  and  Gruev  [19]  (Fig.  30)  would  be  ideal  after  better 
characterization  algorithm  performance  is  achieved  using  the  functional  form  of  r.  Using 
the  Gruev  polarimeter  would  also  allow  for  field  testing  to  ensure  this  algorithm  continues 
to  work  outside  the  lab  environment. 
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